As we begin to construct models, let’s think about the building blocks of mathematical relationships: functions and their shapes.
As a brief review, we will go over the basic function types we will see in this course. Some are trivial, and some are complicated.
The linear function is given by \[ f(x) = ax + b \]
A quadratic function is given by \[ f(x) = bx^2 + cx + d \] A cubic function is given by \[ f(x) = ax^3 + bx^2 + cx + d \]
The exponential function is the unique function \(f(x)\) such that \(f'(x)=f(x)\) for all \(x\). When we write \(e\), it means the “natural number” \(e\approx 2.7182\ldots\). \[ f(x) = a e^{bx} \] The exponential function may also be written as \(f(x) = a \exp(bx)\). When a and b are positive, \(f(x) = a\exp(bx)\) grows quickly. Exponential growth happens when a quantity increases at a rate proportional to how big it already is. More on this later.
The inverse function of \(\exp(\cdot)\) is the natural logarithm \[ f(x) = \log(x) \] That is, if \(x\) is any real number and \(y=\exp(x)\), then \(x = \log(y)\)
Logistic growth happens when initially a quantity grows exponentially but then saturates, or plateaus, as its argument gets large. \[ f(x) = \frac{e^{ax}}{1+e^{ax}} = \frac{1}{1+e^{-ax}} \] The logistic transform takes a real-valued \(x\) and rescales it to be between 0 and 1.
A useful class of sinusoidal functions is \[ f(x) = a + b\sin(c(x-d)) \] where \(a\) is the offset, \(b\) is the amplitude, \(c\) is the rate or frequency, \(d\) is the phase shift. In this course, the argument to all sinusoidal functions will be in radians.
Most models we will see in this class express their dynamics in terms of states of being.
Examples: HIV+, in treatment, using drugs, etc.
For us, a state might be an attribute of an individual person/agent, or the proportion of a population in that state, or simply the state of a single system.
The dynamics of a model describes transitions between states.
A compartmental diagram is a conceptual tool that can be helpful for visualizing model structure and dynamics.
But as we saw above, diagrams may show only part of the story.
Consider a model of health for a single individual, who can be either healthy (\(H\)) or sick (\(S\)). The individual may transition between these states, possibly with different rates.
Suppose \(H=1\) means that the individual is healthy, and \(S=1\) means the individual is sick. Because these are the only possibilities, \(H=1-S\). We say that \(H\) and \(S\) are discrete indicator variables for the individual whose state we are modeling.
Consider a model of health for a population of individuals, each of whom can be either healthy (\(H\)) or sick (\(S\)). The individuals in the population may transition between these states, possibly with different rates.
Let \(H\) be the proportion of individuals who are healthy, and \(S\) the proportion who are sick. Because these population proportions must sum to one, \(H+S=1\).
The compartments are the same, but the level of focus of the model is different: individual versus population. Now \(H\) and \(S\) are continuous instead of discrete.
Discrete-time models index compartmnets by chunks of time: seconds, minutes, days, months, quarters, trimesters, years, follow-up visits, etc. Example: At time points \(t=1,2,3,\ldots\), an individual in \(H\) transitions to \(S\) with probability \(p\), and stays in \(H\) with probability \(1-p\). Likewise, an individual in \(S\) transitions to \(H\) with probability \(q\). This is a discrete-time stochastic model for the individual health state.
Continuous-time models index compartments continuously in arbitrarily small time units. Transitions \(H\to S\) happen in continuous time with rate \(\lambda\), and \(S\to H\) with rate \(\mu\).
Deterministic models always have the same output, given the same inditial conditions.
Stochastic models may have different output, given the same inditial conditions.
Whether stochasticity/randomness matters depends on your goals. If there is heterogeneity in the process you want to model, representing this with randomness in individual attributes may make sense.
One way to construct the dynamics of a compartmental model is to specify the rates of transition between model compartments.
What is a rate? Recall that the derivative of a differentiable function \(f(t)\) at \(t\) is defined as
\[ f'(t) = \lim_{h\to 0} \frac{f(t+h) - f(t)}{h} \]
Intuitively this means that the derivative of a function \(f(t)\) at the point \(t\) is a measure of how quickly the function value changes near \(t\). Equivalently it is the slope of the tangent line of \(f(t)\) at \(t\).
It is often easier to think about the dynamics of a process in terms of its rates of change. It is very common that the solution to a system of differential equations is impossible to write down analytically, but very easy to characterize in terms of its rates of change.
So, from a practical perspective, we can often just write the system in terms of its rates of change, and then use a computer to find the “solution” for the dynamics of the system.
In general, the study of analytic solutions of ODEs is a domain of mathematics. We are going to largely avoid analytical solutions and treat these systems as descriptions that can be solved using a computer. We will see several examples next.
How does compartment occupancy change?
The model is parameterized by its transition rates \(\lambda\) and \(\mu\). What does this mean? Let’s think about it in terms of probabilities
Suppose the current state at time \(t\) is \(H\), so \(H(t)=1\). In a small time \(dt\), the probability of transitioning to \(S\) by \(t+dt\) is \(\lambda \times dt\).
Likewise if the current state at \(t\) is \(S\), so \(S(t)=1\), the probability of transitioning to \(H\) is \(\mu dt\).
This implies that the expected value of \(H\) at \(t+dt\)
\[ \begin{aligned} H(t+dt) &\approx \Pr(\text{stay in }H) H(t) + \Pr(\text{move from $S$ to }H) S(t) \\ &= (1-\lambda dt) H(t) + \mu dt S(t) \end{aligned} \]
Rearranging and taking the limit \(dt \to 0\), we obtain the differential equation
\[ \frac{dH}{dt} = -\lambda H(t) + \mu (1-H(t)) \]
We obtain the differential equation system, which has a population interpretation:
\[ \frac{dH}{dt} = -\lambda H(t) + \mu S(t) \] \[ \frac{dS}{dt} = -\mu S(t) + \lambda H(t) \]
This is the rate of change of the population proportion that is healthy: in a small increment of time \(dt\), a fraction \(\lambda dt\) become sick. However, a fraction \(\mu dt\) of sick people become healthy.
Let \(y(t)\) be the number of bacteria in a dish at time \(t\). Suppose at time \(t=0\), there are \(y(0) = y_0\) bacterium. Each bacterium is immortal, and divides into two with rate \(\lambda>0\). The number of bacteria at time \(t\) is governed by the differential equation \[ \frac{dy}{dt} = \lambda y(t) \] with initial condition \(y(0)=y_0\). The solution is \[ y(t) = y_0 e^{\lambda t} \] You can verify this by computing \(dy/dt\). The dynamics of this system obey “exponential growth”. The speed of this growth is governed by \(\lambda\).
spend some time on this! computing marginal and conditional probabilities on decision trees.
Mechanistic models come with stories connecting physical interactions in real life with the model specification. Make sure your story makes sense.